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Abstract. We study conserved models of crystal growth in one dimension [dtz{x, t) = 
—dxj{x,t)\ which are linearly unstable and develop a mound structure whose typical 
size L increases in time (L ~ t"). If the local slope (m = dxz) increases indefinitely, n 
depends on the exponent 7 characterizing the large m behaviour of the surface current 
j (j - l/Imp): n = 1/4 for 1 < 7 < 3 and n = (1 + 7)/(l + 57) for 7 > 3. 



The conserved dynamics of a solid surface growing under the action of an external 
flux of particles is described by the continuum equation 

dtz{x, t) = -dj{x, t) + 5F{x, t) , (1) 

where z{x,t) is the local height of the surface in a comoving frame (so that the average 
value z is set to 0) and 5F{x, t) is the shot noise. 

Thermodynamic and kinetic mechanisms contribute to j and its actual expression 
depends on the details of the growth process. Here we are interested in the growth of 
a high-symmetry surface by Molecular Beam Epitaxy (MBE), where the instability has 
a purely kinetic origin: the reduced interlayer diffusion Nonetheless, our treatment 
will be as general as possible. 

A wide class of models is described by the current 

i = Km"{x) + JesM , (2) 

where m = d^z is the local slope. The first term generally describes a thermally activated 
relaxation of the surface, but kinetic mechanisms can also contribute to if 0. 

The second term is responsible for the instability and its origin is an asymmetry 
in the sticking process of an adatom to a step (Ehrhch-Schwoebel (ES) effect): sticking 
from the upper terrace is hindered and this implies an up-hill current 0] which is 
called Ehrlich-Schwoebel current (jes)- Also other (generally stabilizing) processes can 
contribute to Jes and this explains the different expressions Jes may take 0]. 
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Whatever these processes are, Jes is hnear in m at small slopes (jes ~ '^^) and 
therefore in the early stages of the growth it prevails on the first term {Km") at 
sufficiently large wavelengths. This means that the linear stability of the flat surface 
will be decided by the sign of u, a positive one meaning instability. In fact, in the limit 
m ^ we have 

dtz = -Kdlz - vdlz (3) 

whose solution is z{x,t) = exp{ujqt) cos(gx) with Uq = vq^ — Kq^. An up-hill current 
means that Jes has the same sign as the slope, so z/ is positive and the fiat surface 



is unstable {uq > 0) against modulations of wavevector smaller than q = ^u/K; the 
instability appears after a typical time of order t* ~ {i^<l'^)~^ = K/u"^. 

The later evolution of the surface depends on the nonlinear form of the unstable 
current jEs{rn). By taking the spatial derivative of Eq. we obtain 

dtm = dl{-j) + d,{6F) (4) 

and a parallel with a phase ordering process is easily done, once we remark that the 
current can be obtained by a pseudo free energy JF: 



|(9,.m)2 + \/(m) 



V'{m) = -jEs(m) .(5) 



The instability of the fiat surface (j'^siO) > 0) means that the potential V{m) has 
a maximum in m = {V"{0) < 0). Contiguous regions of increasing and opposite 
slope are formed. The usual phase ordering process is obtained when V{m) has the 
classical double well form: V{m) = —{v/2)m? + (z//4mQ)m^, corresponding to a current 
Jes = i'm{l — m? / rn^) . After the slope has attained a fraction of mo the dynamics enters 
in the nonlinear regime: the wavelength L of the profile increases in time (coarsening 
process) and the slope saturates to the constant values ±mo. The coarsening law is 
known to be logarithmic ||^ {.L(t) ~ Int) in absence of shot noise and a power law 
{L{t) ~ t^/^), in presence of it. 

The aim of the present paper is to analyze the deterministic {6F{x, t) = 0) growth 
process when V{m) has no minima, corresponding to the absence of zeros at finite slopes 
in the current Jes- We will consider the class of currents defined by 

and the corresponding models will be termed a-models. 

Model-1 has been studied numerically by Hunt et al and they found a coarsening 
exponent n ^ 0.22 {L{t) ~ t"). a-Models without noise have been studied analytically 
by Golubovic through scaling arguments and he finds n = 1/4 irrespectively of a. 
Finally, qualitative considerations based on noise effects give n = l/{2/a + 3), i.e. 
n = 1/5 for model-1. 

Our analytical approach is based on the linear stability analysis of the stationary 
configurations j[m{x)] = 0. This way, the finding of the coarsening exponent n passes 
through the determination of the lowest eigenvalue of the operator (— where H is 
the Hamiltonian corresponding to a particle in a periodic potential ^ . 
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Before proceeding we render adimensional the growth equation by rescahng x with 
1/q, t with t* and z with 1/qi: 

m 

dtz = —dxj , j = m" + — r— . (7) 

Stationary configurations are the solutions of the differential equation j[m{x)] = 
m" + j-Es{fn) = 0. Therefore they correspond to the periodic orbits of a particle in 
the potential —V{m) = — [l/2(a — 1)](1 + m^)-*^"" for a > 1 and in the potential 
—V{m) = (1/2) ln(l +771^) for a = 1. In the former case the potential is upper bounded 
and the solution corresponding to the boundary conditions m —>■ ±oo when x ±oo 
does exist, while it does not for a = 1 because the corresponding energy would be 
infinite. Stationary solutions may be labelled with their period, i.e. the wavelength L: 
rriL^x). 

Let us now perform a linear stability analysis around these stationary and periodic 
solutions: m{x,t) = rriL^x) + ip^x.t). It is easily found that 

dti^ = dl[-r{x,t) + UL{x)i^] , (8) 

where Ul{x) = — jgg(mi(x)). By putting ipi^x^t) = 0(a;) exp(— et) we obtain 

i-dl) [-<f)"{x) + UL{x)<j)] = D,H(t>{x) = e0 . (9) 

Negative eigenvalues mean that mx(x) is linearly unstable and this induces the 
coarsening process; moreover, e(L) 0^ when L oo. The dependence of the ground 
state (GS) energy on the distance L determines the time scale of the coarsening process: 
t ~ l/|e(L)|. For the moment we will assume = 1, i.e. we will consider the 
nonconserved model: dtm = —6 J-'/ 6m. 

First of all we observe that in the limit of large L the energy shift e(L) for the 
periodic potential is equal (up to a numerical factor) to the shift for a single couple of 



potential wells [|T0]. The solution of the problem is given [|TT| in terms of 0o and 0i, 
respectively the ground state for the single well f/i(x), centered in a; = L, and for the 
double well U2{x), centered in a; = ±L. In fact the Schroedinger equations are: 

-0o + f/i0o = O (a) - 0'/ + = (b) (10) 

and by evaluating the quantity /q°° dx[0i x (p!0| a) — (pQ x (p!0^ )] = 0, we obtain 

0i(O)0[,(O) = -e / dxM^)Mx) , (11) 
Jo 

where we have made use of ?7i = f/2 for a; > 0. 

Before proceeding we must determine the asymptotic expressions of 4>o{x) and 
01 (x). The potential U{x) = — jEs(^) is given, for a-models, by 

= ^ . (12) 

^ ^ (1 + m2)"+i ^ ^ 

The asymptotic behaviour of the single-mound profile is obtained by integrating 
the equation m"{x) + jEs('^) = and taking the limit x 00 : 

(Iff) 1 1 

(l/2)(m')^-r(m) = ^ 1 I i.-i • (13) 

dx \Ja — 1 \mf ^ 



Letter to the Editor 



4 



The result m"(x) ~ {a/y/a — l)x, when inserted in (|T2|) gives: 

^ P"- !)("-!) J. , iL , (14) 

with a increasing between a = (for a = 1) and a = 2 (for a = oo). 

The solution of the Schroedinger equation (|To|a) for Ui{x) ~ a/{x — L)^ gives a 
power-law decaying wavefunction {4>o{x) ~ |a; — i^|~^), with an exponent /3 = (1 — 1/a). 

If a < 2 then /3 < 1/2 and therefore the GS 0o(2;) of the single well is not a bound 
state, since da;0o(^) = O'^ the other hand, for a > 2 (po{x) is a bound state and 
01 (x) can be approximated with the expression 0i(x) = [0o(2;) + 0o(^2;)]/\/2. This 
way, from (|Tl]) we easily obtain the relation 

e ~ -20o(O)0o(O) ^ -L-(2/3+i) [a > 2 and = 1] . (15) 

If a < 2, we can put 4>i{x) = [(Pq{x) + 4>o{—x)]/\^ where 0o is a generalization of 
00 to a negative eigenvalue: — 0q(x) + (a/a;^)0o(a;) = e0o(a;). In fact, even if 0o is not a 
bound state, 0i is bounded, because the GS energy e is strictly lower than [/2(±oo) = 0. 
The previous expression for 0i may be used even if 0o itself is bounded (i.e. for a > 2) 
and the result for the coarsening exponent does not change. 

The asymptotic expression for 0o is 0o(2^) = \/^K^{^J\e\x) where is the modified 
Bessel function of order /i = /? — (1/2). The function 0o decays as a power-law 
{4>o{x) ~ lel"^/^"^/^^"^) if a/x^ ^ |e| and exponentially (0o(x) ~ lel""*^/^ exp(— y^x)) 
in the opposite limit, a/x^ <^ |e|. Eq. (p!T| ) now writes 

e / dx0o(x)0o(x) = -20o(O)0[,(O) [a < 2 and D., = 1] , (16) 

JO 

where 0o(x) depends on e. Let us remark that the integral / on the left hand side does 
converge even if 0o is not a bound state. 

The evaluation of the two sides of Eq. (^) is a bit lengthy and we report here 
the result only: |e|ln(l/|e|) ~ l/L^ if a = 2 and |e| ~ l/L^ if 1 < a < 2. In 
Fig. we compare the analytical results for the exponent characterizing the energy 
shift |e(L)| ~ L"^/" with those obtained through its direct numerical evaluation [1^ 
and the agreement is very good. 

Therefore, for the nonconserved model we can conlude that: 

[nonconserved] n = - (1 < a < 2) and n = — {a > 2), (17) 

2 3 — 2/a 

with a logarithmic correction for a = 2 (L ~ (t/ Int)^/^). 

The reason why the coarsening exponent n keeps constant for a < 2 is the following: 
if a > 2 the single-well wavefunction is a bound state, the integral J is a constant 
while the 'superposition' between 0o(x) and 0o(— x) (that is to say the right-hand side 
of (|1^)) decreases at increasing a, which implies a decreasing n. Conversely, when 
a < 2 the integral / becomes a— dependent and decreases with a: these dependence 
counterbalances the reduction of the right-hand side of (0). 

For the conserved growth model, Dr^ = —d^ and Eq. (|ll|) must be replaced by a 
more complicated expression. It has not been possible to carry out a rigorous calculation 
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because [0i x DxH(j)Q — (pQ x D^Hcf)]] is no more integrable. Nonetheless, there are 
strong indications that the right-hand sides of (p^JT6| ) acquire a factor the origin 
of this scahng factor is that (t)o{x) has a power-hke behaviour (and therefore derivation 
corresponds to divide by x) and also that U{x) ~ x~'^ . Furthermore, since we need the 
single well wavefunction, corresponding to a zero energy, a solution of the Schroedinger 
equation H(j){x) = is also solution of Dj.H<t>{x) = 0. 

As a consequence of such factor, the coarsening exponent for the conserved case is 
easily obtained from the nonconserved one: (1/n) [(!/'"') + 2]- Therefore: 

[conserved] n = - (1 < a < 2) and n = — (a > 2) . (18) 

^ ' 4 ^ ~ ^ 5-2/a ^ / V / 



In order to check numerically the validity of the results reported in Eg .([181) 
and therefore the dependence of the coarsening exponent n on the parameter a, 
detailed numerical simulations have been performed. In particular, we have numerically 
integrated equation (|^ by employing a pseudo-spectral time splitting code . 



The values of L(t), whose log-log plot gives the exponent n, are evaluated 
through the power spectrum (PS) of z{x,t): the weighted average of the wavevectors 
corresponding to the most relevant components of the PS is 2n/L{t). A different method 
using the spatial correlation function gives consistent results. In Fig. ^, the numerical 
findings for n{a) by direct integration of Eq. (|^) are shown together with the theoretical 
expression (|18|) and a good agreement is found. 

In conclusion we have found the analytic expression for the coarsening exponents 
n{a), both for the nonconserved model, Eq. (|1^) and for the conserved one (growth 
model), Eq. (|TB|). Coarsening varies with a and it is not logarithmic (i.e. n = 0) even 
for a = oo. 
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Figure 1. Analytical (full line) and numerical (crosses) values for the exponent 1/n 
governing the asymptotic energy shift |e2| ~ l/L^^" (nonconserved model). 



Figure 2. Coarsening exponent n for the conserved model. In the inset we enlarge 
the small a region. Full line is the analytical result (Eq. |l^ . Points are the exponents 
found integrating numerically Eq. ^ for a system size M = 1024 (spatial resolution 
Ax = 0.25) and a total time 400, 000 < T < 1, 600, 000 (time step At = 0.05). A few 
tests have also been done with a smaller time step {At = 0.025) and longer chains 
(M = 2048-4096), obtaining consistent results. Bars indicate the numerical fit errors. 
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